Metastability of persistent currents in trapped gases of atoms 
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We examine the conditions that give rise to metastable, persistent currents in a trapped Bose- 
Einstein condensate. A necessary condition for the stabihty of persistent currents is that the trapping 
potential is not a monotonically increasing function of the distance from the trap center. Persistent 
currents also require that the interatomic interactions are sufficiently strong and repulsive. Finally, 
any off-center vortex state is shown to be unstable, while a driven gas shows hysteresis. 
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I. INTRODUCTION 

Vapors of ultracold atoms provide an ideal laboratory 
for testing the superfluid properties of quantum systems, 
including - among many other things - persistent cur- 
rents. These phenomena are not only interesting from a 
theoretical point of view, but may also have very impor- 
tant technological applications [ij. 

As opposed to the traditional studies on liquid Helium 
that is confined in a box, gases of atoms may be confined 
in potentials of various form, such as harmonic, anhar- 
monic, etc. In addition, the coupling constant between 
the atoms is tunable, and it is even possible to change 
its sign, allowing us to study both repulsive, as well as 
attractive effective interactions between the atoms. 

As we demonstrate below, the two above features allow 
one to tune the stability of persistent currents in atomic 
systems by appropriately choosing the functional form of 
the trapping potential, the size/strength of the scattering 
length that describes the low-energy two-body collisions 
between the atoms, or the atom density/atom number. 

Numerous theoretical and experimental studies have 
examined the rotational properties of clouds of trapped 
atoms, the vortex nucleation and stability, and more gen- 
erally the superfluid properties of these gases. For a re- 
view of the extensive work on all these problems, we refer 
to the two articles by Leggett, on superfluidity in general 
0, and in trapped gases [3]. 

In the present study, we examine the metastability 
of a current-carrying state with one unit of circulation, 
or equivalently, the energetic stability of a vortex state 
that is located at the center of the cloud. Our results 
show that in a two-dimensional trap, at least in the 
Thomas-Fermi limit of strong interactions, the vortex 
state, viewed as a particle, feels a force that is propor- 
tional to the gradient — Vno(/9), where no(p) is the single- 
particle density distribution of the non-rotating cloud. 
Here p is the usual radial variable in cylindrical polar co- 
ordinates. Therefore, we conclude that a necessary con- 
dition for metastability is that 7io(p) is not a monotoni- 
cally decreasing function. We give an analytic derivation 
of this result in Sec. VI. 

For weak interactions metastability is not possible. 
The dispersion relation, i.e., the total energy of the cloud 



£ for a given value of the angular momentum per atom 
Ih, is dominated by the single-particle energy, and there- 
fore local minima cannot appear in £{l). In the Thomas- 
Fermi limit of strong and repulsive interactions, metasta- 
bility is possible, since the interaction energy is much 
larger than the single-particle energy. In this limit, the 
density nQ{p) of the non-rotating cloud is a mirror image 
of the trapping potential V{p). Therefore, we conclude 
that in the Thomas-Fermi limit a necessary condition for 
metastability is also that V{p) is not a monotonically 
increasing function of p. 

In our discussion so far we have distinguished be- 
tween ?T.o(p) and V{p). However, since metastability 
requires that the interaction energy is at least compa- 
rable to the single-particle energy, and in this case a 
monotonically-increasing V{p) implies a monotonically- 
decreasing uq (p) , a necessary condition for metastability 
is also that V{p) is not a monotonically increasing func- 
tion of p. 

In what follows we first present our model and our 
method in Sec. II. Section HI describes the results we 
have derived within the mean-field approximation, and 
Sec. IV the results of numerical diagonalization. Section 
V examines the case of attractive interactions. In Sec. VI 
we give some analytic arguments which support our nu- 
merical results. Finally, in Sec. VII we consider a driven 
gas, and conclude that independently of the form of the 
trapping potential, any single off-center vortex state is 
unstable, while the gas shows hysteresis. 



II. MODEL AND METHOD 

In general one has to solve a three-dimensional prob- 
lem. However, the motion along the axis of rotation 
introduces an extra degree of freedom, which may give 
rise to vortex bending. To get rid of this complication, 
we thus consider a highly-oblate trap, with a very tight, 
harmonic confinement along the axis of rotation - taken 
to be the z axis - and some axially-symmetric trapping 
potential V{p) perpendicular to the axis of rotation. 
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Kxt(r) = Vip) + -Mlu 



(1) 
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We assume that the corresponding quantum of energy 
for motion along the z axis is much larger than any 
other energy scale in the problem. The motion of the 
atoms is thus quasi-two-dimensional, since their motion 
along the z axis is frozen. 

The energy of the gas a.t I — 1, £{l = 1) is higher than 
£{l = 0). The interesting question is whether there is a 
barrier in the dispersion relation that separates the 
state with I = 1 from the state with 1 = 0. In the presence 
of such a barrier, the current-carrying state is metastable, 
and relatively weak perturbations that do not conserve 
the energy and the angular momentum cannot destabilize 
it. 

We mentioned above that our study examines the en- 
ergetic stability of a vortex state that is located at the 
center of the cloud. To understand the connection be- 
tween this kind of stability criterion, and the calculation 
of one should remember that there is an one-to- 

one correspondence between the angular momentum per 
atom Ih, and the distance of the vortex state from 
the center of the trap. An expression for l{Ry) is de- 
rived in Sec. VI. The function £{l) may thus be viewed 
also as a function £{Rv), and since l{Rv = 0) = 1, a 
local minimum in £{l) for / « 1, implies that the vortex 
state is energetically stable under (at least) infinitesimal 
displacements of the vortex from the center of the trap. 

To calculate the dispersion relation, we use both the 
mean-field Gross-Pitaevskii approximation, as well as nu- 
merical diagonalization of the Hamiltonian. The calcu- 
lated energies from the two methods agree to leading or- 
der in iV [J|, i.e., they differ due to finite-size corrections 
in the small systems that we consider in numerical diag- 
onalization. 



III. MEAN-FIELD APPROXIMATION 

Starting with the mean-field, Gross-Pitaevskii approxi- 
mation, the form of the trapping potential (which is tight 
along the z axis) allows us to use the following ansatz for 
the order parameter, 

■i>{x,y,z) = <^{x,y)M^'), (2) 

where (f>o{z) = e~^^^^'^°'''> /{TralY^"^ is the ground state of 
the harmonic oscillator along the z axis, with Uz being 
the oscillator length along this axis. The corresponding 
Gross-Pitaevskii equation for the two-dimensional order 
parameter ^(x,y) is 

where ^ is the chemical potential, and g = 
NUoJ\Mz)\^dz = Vs^h^Na/Maz. Here N is the 
number of atoms, and Uq = ^Trfi^a/M, where a is the 
s-wave scattering length for elastic atom-atom collisions, 
and M is the atom mass. We define for convenience the 
dimensionless quantity g = Mg/fi^ = \/8TTNa/az. This 



parameter may easily get at least as large as 10^, if, for 
example, N = 10^, a = 50 A, and ^ 10 fim. 

The functional whose minimum value we are looking 
for, is 

= -|g / ^*iyl + yl)<^dxdy + 

J ^*V{x,y)^dxdy + ^ J l^l^dxdy, (4) 

under the extra constraint of a fixed normalization 
J |<I>P dxdy = 1. Typically one minimizes the above ex- 
pression at a fixed rotational frequency of the trap ft, 
by introducing the new function ^rot = £ — l^, where 
I = i J ^*{ydx — xdy)^ dxdy. If the second derivative of 
the (minimized) function £ with respect to / is positive, 
then one may get the minimized value of £{Q), as well as 
at fixed ri, and finally derive the dispersion relation 

£{l). 

However, if the second derivative of the (minimized) 
function £ with respect to I is negative, within a given 
range of I, li < I < I2, the only minima occur at the end 
points I = h and I = 12- To overcome this difficulty, one 
may minimize the following expression instead [5j, 

E{^\<^)=£{^*,^) + ^hu{l-h)\ (5) 

again under the constraint of a fixed normalization 
J <^*^dxdy = 1. Here C and /q are real, dimensionless, 
and positive constants. The constant C has to be chosen 
sufficiently large to ensure that the last term in Eq. ([5]) is 
the dominant one, which then also gives the correspond- 
ing minimized function E{1) a positive curvature, namely 
£"{l) + Chcu > 0, for all values of the angular momentum 
h < I < I2, where £"{l) denotes the second derivative 
with respect to Under the condition of a positive cur- 
vature of E{1), one may minimize E(l) at any Zmin that 
satisfies the equation 

Un = /O - £'{lrain)/iChL0). (6) 

Using the method of imaginary-time propagation 
we minimize the expression of Eq. ^ for three different 
trapping potentials, namely harmonic, anharmonic, and 
"Mexican-hat-hke" , that are plotted in Fig. 1. In Fig. 2 
we consider the harmonic trapping potential of the form 

Vip) = \mu:^p\ (7) 

and plot the dispersion relation £(?), in units of Tiuj, for 
four values of the parameter g = 10, 100, 500, and 1000. 
In this figure we have shifted the zero of the energy, so 
that for all curves £{l = 0) = 0. The actual shifts of the 
energies are: for g = 10, £{l = 0) « 1.59?ia;, for g = 100, 
£[l = 0) «:! 3.95fi.w, for g = 500, £[l = 0) « ^.blhw and 
for g = 1000, £{l = 0) = ll.Qlhw. 

As one can see in Fig. 2, the curvature of £{l) is nega- 
tive for all values of g. Although the slope of this curve 
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FIG. 1: The three forms of the trapping potential V{p) that 
we use in the present study, namely harmonic, Eq. ((Tjl, an- 
harmonic, Eq. and "Mexican hat", Eq. The energy is 
measured in units of Hlj, and the length in units of a±. 



at I = 1~ decreases with increasing (7, our calculations 
show that it never becomes zero, and no metastable min- 
imum forms around 1 = 1. Therefore, independently of 
the value of the coupling, a harmonic trapping potential 
cannot support persistent currents. 

To examine the effect of the functional form of V{p), 
we also consider a more flat potential around the origin, 
which has the form 



V = ^huj{p/a^r, 



(8) 



where a_L = (?i/Mcj)^/^ is the oscillator length along the 
plane of motion of the atoms, for a = 32, and g = 30, 100, 
and 1000. Actually, this high value of a corresponds to 
a trap which is almost like a hard-wall potential, with a 
radius a±. 

Again, we plot £{l) in Fig. 3 shifting the zero of the 
energy in each case, so that £{l = 0) = 0. These shifts 
are: for g = 30, £{l = 0) « 7Alhuj, for g = 100, £{L = 
0) w 67.15huj, and for g = 1000, £{L = 0) = 124:.72?uj. 

This almost hard-wall potential does not allow the gas 
to expand radially, as in the case of harmonic confine- 
ment. The slope of the dispersion relation £{l) for Z ^ 
increases with increasing g; in contrast for a harmonic 
trapping potential, the slope decreases. This observation 
probably implies that there is a critical power-law depen- 
dence of the trapping potential, for which this slope does 
not depend on g. Also, the slope of £{l) for / 1~ tends 
to zero for high values of g, due to the almost homoge- 
neous density of the cloud close to the center of the trap. 
Apart from these observations, £{l) has the same qualita- 
tive features as in the case of harmonic trapping, with no 
metastable minimum forming around 1 = 1. Persistent 
currents are not possible in this case, either. More gen- 
erally, our calculations imply that this is also true even 
for a hard- wall potential, a ^ 00. 
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FIG. 2: Top: The dispersion relation £{V) in units of fiu), for 
a harmonic potential, and for g = 10, 100, 500, and 1000, cal- 
culated within the mean-field approximation. Bottom: E{V) 
that results from numerical diagonalization of the Hamilto- 
nian (dots), for A'^ = 4 atoms and < L < 6 units of angular 
momentum, as compared to the one calculated within mean 
field, for g = 2. 



We turn now to another functional form of V{p)^ which 
does not increase monotonically with p, i.e., we consider 
a "Mexican-hat" shape. 
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V{p) = -Mlu^p^ + Voe 



2 / 2 



(9) 



with Vq > hLu/2. Such a potential can be realized exper- 
imentally if one applies a Gaussian, narrow laser beam 
on top of the ordinary harmonic trap Q ■ 

Figure 4 shows the corresponding dispersion relation, 
again after shifting the zero of the energy at Z = 0, for 
Vo/hcj = 4, and 5 = 2, 20 and 100. The shifts are: for 
g = 2, £{l = 0) « 2.44^^, for g = 20, £{l = 0) « 2.94^^;, 
and for g = 100, £[l = 0) = A.bQhuj. In this case £{l) 
does develop a local minimum around ^ = 1 for sufh- 
ciently strong interactions, g ~ 20. Persistent currents 
are possible in this case. In Ref. Q we have shown that 
in a strictly one-dimensional trap with periodic bound- 
ary conditions, a similar picture emerges for £{l) as the 
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Angular momentum per atom 

FIG. 3: The dispersion relation £(l) in units of huj, for an 
anharmonic potential of the form of Eq. (|8]), with a — 32, and 
for g = 30, 100, and 1000, calculated within the mean-field 
approximation. 

strength of the interaction increases. This is not a sur- 
prising result, since, as shown in Fig. 2, the trapping po- 
tential is - roughly speaking - toroidal-like in this case. 

IV. NUMERICAL DIAGONALIZATION 

We have also diagonalized numerically the many-body 
Hamiltonian, for small systems, of N atoms and L units 
of angular momentum. Again we assume tight confine- 
ment along the axis of rotation and essentially consider 
two-dimensional motion. We use the eigenstates of the 
harmonic oscillator in two dimensions to build the Fock 
states, which are eigenstates of the operators N and L, 
and then diagonalize the resulting Hamiltonian matrix. 

Clearly one has to truncate, since the dimensionality 
of the matrix becomes very large if one considers highly- 
excited states. In the present calculation for each L and 
N, we include in the Hilbert space as many single-particle 
states as are necessary, so that the lowest eigenvalue of 
the Hamiltonian has converged up to the third decimal 
point. 

The calculated energies that we obtain agree with 
those of the mean-field approximation, apart from finite- 
size corrections {5]. The dots in the lower plots of Figs. 2, 
and 4 show the corresponding dispersion relation (i.e., the 
lowest eigenenergy divided by N), as function of the an- 
gular momentum per atom L/N, for < L < 4, = 4, 
and 5 = 2. 

V. ABSENCE OF METASTABILITY FOR 
ATTRACTIVE INTERACTIONS 

Another interesting question is the possibility of 
metastable currents in the case of an effective attractive 
interaction between the atoms. Our calculations show 




Angular momentum per atom 

FIG. 4: Top: The dispersion relation £{l) in units of hui, for a 
"Mexican-hat" potential of the form of Eq. Q, with Vo/hoj = 
4, and for g = 2, 20, and 100, calculated within the mean-field 
approximation. Bottom: £{l) that results from numerical 
diagonalization of the Hamiltonian (dots), for A'^ = 4 atoms 
and < L < 6 units of angular momentum, as compared to 
the one calculated within mean field, for g = 2. 



that metastability is not possible in this case, in agree- 
ment with the arguments of Leggett 3 . 

In two interesting studies, Wilkin, Gunn, and Smith 
[3], as well as Mottelson [ifl] have shown that in a har- 
monic trap, for an effective interaction that is very weak 
and attractive, within the lowest Landau level approxi- 
mation, the angular momentum is carried by the center of 
mass motion ^11]. The reason is that in a harmonic trap, 
the center of mass and relative motions decouple, and 
as a result the least energetically-expensive way to give 
angular momentum to the gas, is via excitation of the 
center of mass. The relative distance between the atoms 
is thus unaffected, and the interaction energy does not 
depend on the angular momentum. 

More generally, when the effective interaction is at- 
tractive, but not infinitesimally weak, the lowest, non- 
rotating state is already an interesting and non-trivial 
question. The absolute minimum of the energy of the 
gas for negative g corresponds always to a collapsed 
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FIG. 5: The dispersion relation f (i) in units of Ivjj, calculated 
within the mean-field approximation, for a "Mexican-hat" po- 
tential of the form of Eq. ((2)), with V^lhix) = 4, and for g = —2 
(higher curve), and g = —3 (lower curve). 



state. However, before the cloud collapses, for sufficiently 
small values of the parameter \g\, the cloud may be in a 
metastable state. If the trapping potential V{p) increases 
monotonically with p, the single-particle density of the 
(non-rotating) cloud is axially-symmetric, but it shrinks 
towards the origin of the trap, because of the attractive 
interaction. If the trapping potential does not increase 
monotonically, the gas may have a homogeneous density 
along the axially symmetric minimum of the potential, 
or it may form a symmetry-breaking localized blob. 

Depending on the lowest state in the absence of rota- 
tion, the current-carrying states may involve vortex exci- 
tation, or center of mass excitation. Our results indicate 
that, independently of the form of the trapping poten- 
tial, and of the single-particle density distribution, the 
dispersion relation is a smooth function of with a non- 
negative curvature. 

Figure 5 shows a specific example, where we calculate 
£{l) for a trapping potential that has a "Mexican-hat" 
shape, with Vq/Tilo = 4. The density is axially symmet- 
ric around the minimum of V{p) for g — —2, and it is 
localized for g = —3, having formed a soliton-like blob, 
as in the case of the one-dimensional problem [l^, [isj . 
As g becomes even smaller, the cloud collapses. In both 
cases (7 = — 2, and = — 3, the dispersion relation has a 
positive curvature. In particular when g — —3, £{l) cx P, 
since the motion corresponds to solid-body rotation. 



possible. We give here a simple derivation of this result. 

The dispersion relation £{l) consists of the single- 
particle energy and of the interaction energy, £ — {Hsp) + 
{Vint}- The contribution of the single-particle part of 
the Hamiltonian Hsp is a smooth, and increasing func- 
tion of /, with d{Hsp)/dl > 0. We thus focus on the 
interaction energy (Mnt), i-e., on the quantity (T4nt) = 
{g/2)J\<^\Uxdy. 

We consider the derivative 



dl 
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dl 



■ dxdy. 



(10) 



Clearly, the sign of the above derivative is determined by 
the sign of The density that enters in the above 

expression is the local density of the cloud 7T,(p, (f) = 
In the limit of strong interactions, g 3> 1, the coherence 
length ^ is much smaller than the radius of the cloud 
i?o, with ^/i?o ^ .9^^^^- The distortion of the cloud 
(due to the vortex) is localized around a region of radius 
« ^. For this reason, in the calculation of the integral 
in Eq. PTI| one may use the density profile neglecting the 
density depression due to the vortex, and also the weak 
dependence of the density of the gas on the polar angle 
<f>. Denoting the corresponding density as n{p) [which is 
qualitatively the same as no(p)]. 



-^^g n{p)^dxdy. 



Then, using the chain rule 



dn dRy 
_ = (p.Vn) — 



(11) 



(12) 



where Ry is the distance of the vortex from the center 
of the cloud, we observe that when n{p) decreases mono- 
tonically, the inner product p ■ Vno is always negative. 

Turning to the derivative dRy/dl, this is negative. The 
function Ry{l) is monotonically-decreasing, with Ry — 
for I = 1, and Ry = oo, when I = 0. Generalizing the 
result of Ref. [l3| that is given for a homogeneous system, 
we derive an approximate expression for l{Ry). More 
specifically the angular momentum of the gas around the 
center of the trap is (in units of fi,). 



— / AIpv n{p) pdpdcj). 



(13) 



However, the circulation 



VI. ANALYTICAL ARGUMENTS 

The results of Sees. HI, IV and V show that metasta- 
bility of superfiow is possible for sufficiently strong and 
repulsive interactions. Furthermore, when the single- 
particle density no {p) of the non-rotating cloud decreases 
monotonically, £{l) is a monotonically-increasing func- 
tion and local minima in the dispersion relation are not 



V • dw = 



h_ 



(14) 



provided that the vortex is inside the area that is de- 
fined by the corresponding line w. Considering a circle 
of radius p, 



(15) 
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for p> Ry and zero for p < From Eqs. (fT3|) and (fTSl) 
we get that, 

POO /"OO 

l{Ry) = 2TT n{p)pdp= n{p)d^p, (16) 

which is the fraction of the atoms that reside outside a 
circle of radius Equation implies that l{Rv) is 
a decreasing function, and therefore dRy/dl < 0. An 
analytic expression for l{Ry) may be derived for weak 
interactions (i.e., within the lowest-Landau level approx- 
imation), which is l{Rv) « 1 — Rl/{2a\) [T3|, when the 
vortex state is close to the center of the trap (for small 
values of i?t,). 

We conclude from Eqs. (fTTj) . (fT2|) and (fT6|) . that when 
no{p) decreases monotonically, d{Vint)/dl > 0, and thus 

which means that indeed metastability is not possible in 
this case. 

As mentioned also earlier, in the Thomas- Fermi limit 
one may go even further with this argument, since no(p) 
is roughly the mirror image of the trapping potential 
V{p). Therefore, in the Thomas-Fermi limit, a necessary 
condition for metastability of superflow is that V{p) does 
not increase monotonically with p. This condition holds 
even for more weak interactions, i.e., in the cross-over 
region, when the single-particle energy is comparable to 
the interaction energy. 

VII. THE BEHAVIOR OF A DRIVEN GAS AND 
SOME GENERAL CONCLUSIONS 



and highly-oblate trap, it has been shown that in the 
Thomas-Fermi limit, n2 ~ (5/3)rii [l^. As we explain 
below, the fact that fii < il2 < ^3 has some important 
implications. 

Experimentally one may either cool down the gas 
to very low temperatures and then rotate, or first ro- 
tate above the condensation temperature, and then cool 
down. In the first case, as one starts rotating the trap (at 
a fixed g), the cloud undergoes a discontinuous transition 
from the non-rotating state to the state with one vortex 
that is located at the center of the trap, at the frequency 
il = ila. In the reverse process, the vortex leaves the 
cloud discontinuously at an — Qi, and therefore, there 
is hysteresis. In addition, for fli < Q < Q2, the vortex 
state is metastable, i.e., it becomes locally stable before 
it becomes thermodynamically stable (at = )• Fi- 
nally, because of the negative curvature of £{l), any single 
off-center vortex state is unstable. 

In the second case, if one first rotates and then cools 
down, a (meta)stable vortex state will appear at the cen- 
ter of the cloud at a rotational frequency = Jli. Again, 
any single off-center vortex state is unstable. 

The physical picture that underlies the above argu- 
ments for the metastability, the hysteresis, and the fact 
that any off-center vortex state is unstable, is that the 
vortex state creates a node at the density of the gas. As 
long as the density no{p) decreases monotonically, this 
node cannot produce any energy barrier [iTj . On 
the other hand, if the density is not a monotonically- 
decreasing function of p, there may be an energy barrier 
in the dispersion relation, provided that the interaction 
is repulsive and sufficiently strong. 



In addition to the question of metastability of super- 
flow, the calculated £{l) for the various forms of V(p) 
also allow us to make some general statements for the 
behavior of a driven gas. Remarkably these statements 
are independent of V{p). 

More specifically, one may define three critical rota- 
tional frequencies, namely the slope fii of £{l) at I — 1^, 
the frequency 1)2 at which h^l2 = £{l = 1) — £{l — 0), 
and the slope fl^ of £{l) at Z = 0. According to our cal- 
culations, in all cases, Qi < Q2 < ^3- In a harmonic. 
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